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Abstract. Systems with long-range interactions, while relaxing towards equilibrium, 
sometimes get trapped in long-lived non-Boltzmann quasistationary states (QSS) which 
have lifetimes that grow algebraically with the system size. Such states have been observed 
in models of globally coupled particles that move under Hamiltonian dynamics either on a 
unit circle or on a unit spherical surface. Here, we address the ubiquity of QSS in long- 
range systems by considering a different dynamical setting. Thus, we consider an anisotropic 
Heisenbcrg model consisting of classical Heiscnbcrg spins with mean-field interactions and 
evolving under classical spin dynamics. Our analysis of the corresponding Vlasov equation 
for time evolution of the phase space distribution shows that in a certain energy interval, 
relaxation of a class of initial states occurs over a timescale which grows algebraically with 
the system size. We support these findings by extensive numerical simulations. This work 
further supports the generality of occurrence of QSS in long-range systems evolving under 
Hamiltonian dynamics. 



PACS numbers: 05.20.-y, 05.70.Ln, 05.40.-a 



1. Introduction 

Long-range interacting systems have attracted considerable interest in recent years [TH6]. 
These systems are characterized by an interparticle potential which in d dimensions decays 
at large separation, r, as l/r a , with a < d. Examples include non-neutral plasmas [7J, dipolar 
ferroelectrics and ferromagnets [Hj, self-gravitating systems [9], two-dimensional geophysical 
vortices [10J, etc. Long-range interactions lead to non-additivity, whereby thermodynamic 
quantities scale superlinearly with the system size. This may result in equilibrium properties 
which are unusual for short-range systems, e.g., a negative microcanonical specific heat 
[TT|[T2]. inequivalence of statistical ensembles [l3l[T3], and many others [T5] . 

As for the dynamical properties, long-range systems often exhibit broken ergodicity 
[HJ[T6] an d intriguingly slow relaxation towards equilibrium [T0 | [T ^ [T7H20] . Such slow 
relaxation has been discussed in the context of self-gravitating systems (see, e.g., [2T||22]). 
and has recently been demonstrated in a model of globally coupled particles (inertial rotors) 
moving on a unit circle and evolving under deterministic Hamilton dynamics. In this so- 
called Hamiltonian mean-field (HMF) model, it has been found that for a wide class of initial 
distributions, the relaxation time diverges with the system size [17]. For some energy interval, 
the relaxation proceeds through intermediate long-lived quasistationary states (QSS). These 
non-Boltzmann states exhibit slow relaxation of thermodynamic observables, and have a 
lifetime which grows algebraically with the system size [HI [23]. As a result, the system 
in the thermodynamic limit never attains the Boltzmann-Gibbs equilibrium but remains 
trapped in the QSS. At other energies, however, relaxation occurs faster on a time which 
grows logarithmically with the system size. Quasistationary states in the HMF model exhibit 
interesting features like anomalous diffusion and non-Gaussian velocity distributions [2H426] . 
Generalization of the model to include anisotropy terms in the energy [23], and to particles 
which are confined to move on a spherical surface rather than on a circle [27] , have also shown 
the existence of QSS. All these models evolve under deterministic Hamiltonian dynamics. 
The robustness of quasistationarity to stochastic dynamical processes has also been analyzed, 
where it is found that QSS exist only as a crossover phenomenon. Under such dynamics, 
these states have a finite relaxation time which is determined by the rate of the stochastic 
process [28H33] . 




In this work, we address the issue of ubiquity of QSS in long-range systems by examining 
a different dynamical model for its existence than the HMF model. To this end, we consider 
an anisotropic Heisenberg model with mean-field interactions. The model comprises globally 
coupled three-component Heisenberg spins evolving under classical spin dynamics. The 
dynamics is thus very different from particle dynamics of previously studied long-range 
systems that have shown the existence of QSS. 

Our model has an equilibrium phase diagram with a continuous transition from a low- 
energy magnetic phase to a high-energy non-magnetic phase. Analysis of the Vlasov equation 
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for evolution of the phase space distribution shows that, as in previously studied models, 
our model exhibits relaxation over times that grow either logarithmically or algebraically 
with the system size. At low energies, a non-magnetic initial state is dynamically unstable. 
It relaxes to the stable magnetically-ordered state on a logarithmic timescale. At higher 
energies, but within the magnetic phase, the non-magnetic state becomes linearly stable, 
and its relaxation to equilibrium occurs on an algebraically growing timescale so that QSS 
are observed. These results, obtained by analyzing the marginal stability of the Vlasov 
equation, are supported by extensive numerical simulations. This analysis yields further 
evidence for the occurrence of QSS in generic long-range systems. 



2. The model 



We start by defining the model of study. It consists of N globally coupled classical 
Heisenberg spins of unit length, Sj = (Si X , Si y , Si Z ), i = 1, 2, . . . , N. In terms of spherical 
polar angles 0, G [0, ir] and <pi G [0, 2ix] for the orientation of the i-th spin, one has 
— sin 6i cos 0j, Siy = sin 0j sin <f>i,S iz = cos 6^. The Hamiltonian of the model is given 

by 

j N N 

^ = -^£ s - s ^£^' (!) 

i,j=l i=l 

where the first term with J > describes a ferromagnetic mean-field like coupling, while 
the last term gives the energy due to a local anisotropy. We take D > such that at 
equilibrium, the energy is lowered by having the magnetization m = (1/N) YliLi S« pointing 
in the xy plane. The coupling constant J is scaled by 1/N to make the energy extensive, in 
accordance with the Kac prescription [33]. However, the system continues to remain non- 
additive in the sense that it cannot be trivially divided into independent macroscopic parts, 
as can be achieved with short-range systems. In this work, we take J = 1 and the Boltzmann 
constant ks = 1- 

The time evolution of the model (CQ) is governed by the set of first-order differential 
equations 

^ = {S,,iJ}; z = l,2,...,iV. (2) 

Here the Poisson bracket {A, B} for two functions of the spins is obtained by noting 
that suitable canonical variables for a classical spin are and S z , so that in our model, 
{A, B} = ^^(dA/dfadB/dSu - dA/dS iz dB /dfa). It may be rewritten as [35] 

r , ^ \— > ^, dA dB , . 

{AB} = £ S ,.-x-. (3) 

i=l 1 1 

Using the above relation, we obtain the equations of motion of the model: 
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2DS iy S lz , (4) 
S ix m z + 2DS ix S iz , (5) 
S iy m x . (6) 

Note that the above equations of motion may also be obtained by considering the 
corresponding quantum equations and taking the limit of infinite spins. 

From Eq. (jHJ), one finds by summing over i, that m z is a constant of motion. The 
motion also conserves the total energy and the length of each spin. From Eqs. (j3J), (jSJ), and 
(jHJ), the time evolution of the variables Oi and (pi are obtained as 

Oi = m x sin (pi — m y cos <pi, (7) 
4>i — m x cot 0, cos (pi + m s cot 9i sin 0j — m z + 2D cos 0j. (8) 

Note that the Hamiltonian of the previously studied models of particles moving either 
on a unit circle [IT] or on a spherical surface [27] may also be expressed in terms of spin 
variables. However, the dynamics of these models is rather different from that of the model 
considered here. Unlike the present model, all Poisson brackets of the spin variables of 
these models are set to zero, and the dynamics is generated by a kinetic energy term which 
is explicitly included in the Hamiltonians and which is absent in the present model. It is 
worthwhile to point out that even in the large D limit, when the ^-component of the spins 
becomes vanishingly small, the dynamics is different from the HMF model with particles 
moving on a unit circle. 



Six SiyTfl z 
Siy Si z Ul x 
Siz Si x TTly 



3. Equilibrium phase diagram 



We now discuss the equilibrium phase diagram of our model. The canonical partition function 
Z is given by 

-/3Nm 2 



Z = / Y\ si* 1 @idOid(pi exp 



13 D cos2 Oi 



(9) 



where j3 is the inverse temperature, and m = (m x + w? y + m 2 z ) 1 / 2 . Since Eq. ([T]) describes a 
mean-field system, it is straightforward to write down expressions for equilibrium quantities 
like the average magnetization or the average energy. As mentioned above, for D > 0, the 
system orders in the xy plane. Choosing the ordering direction as x without loss of generality, 
the average equilibrium magnetization along x, given by (m x ) = (sin 9 cos 0), is 



(m x ) 



f sin 2 9 cos 



,0(m x ) sin cos <fr—/3D cos 2 8 



sin 



gjS K) sin 9 cos (f>—/SD cos 2 9 



(10) 
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Close to the critical point, one may expand the above transcendental equation to leading 
order in (m x ) to get 



e -^ Dcos2e \ =0. nr 



The transition temperature f3 c is obtained by equating the bracketed expression to zero, and 
may be seen to satisfy 

2_ _J_ e~^ D 

& ~ 2f3 c D + v^DErf [y/fcD\ ' (12) 
Here, Erf[x] = J= J* e - * 2 dt is the error function. The critical energy density, e c = D(cos 2 9), 



is 

^(l"£), (13) 

see Fig. (pQ). In model (JTJ), where the phase transition is continuous, the canonical and 
microcanonical ensembles are expected to be equivalent [T5|l36], and Eq. ([13]) therefore also 
yields the microcanonical energy at the transition. 

4. Relaxation towards equilibrium 

To study the dynamical behavior of magnetization, we now analyze the Vlasov equation for 
evolution of the phase space density [TJ. For any mean-field model like ours, this equation 
faithfully describes the iV-particle dynamics for finite time and in the limit N — > oo [3ll37j. 
For the model ([I]), we study the dynamics by examining a single particle, which is moving 
in the two-dimensional phase space of its canonical coordinates, <fi and S z = cos 6, due to 
the mean-field produced by all other particles. Here, and in the following, the particle index 
is suppressed for brevity. Let g(cos9,(j),t) be the probability density in this single-particle 
phase space, such that g(cos 9, <j), t)dcos6d<j) gives the probability to find the particle with 
the z-component of its spin between cos 9 and cos 9 + d cos 9 and the azimuthal angle between 
4> and (p + d(p at time t. In terms of the canonical coordinates, <fi and cos#, the flow in the 
phase space is divergence free. Conservation of probability then implies vanishing of the 
total time derivative of g, that is, 

dg dg ■ dg -dg 

dt at a cos 9 ocp 

where dot represents derivative with respect to time. More conveniently, we define a new 
function, f(9,<p,t) = g(cos9,(f),t), such that / (9, (f>,t) sin 9 d9d<p gives the probability to find 
the particle with its angles between 9 and 9 + d9 and between and + d<p at time t. 
The equation for the evolution of / is straightforwardly obtained from that of g to yield 
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df/dt + 9{df/d9) + j>(df/d<j>) = 0. Using Eqs. © and © to substitute for 9 and gives 
the Vlasov equation for time evolution of f(9, 0, t) as 



df 
dt 



m y cos < 



m x sin - 



9/ 
<9# 



m T cot # cos i 



+ m y cot sin i 



m 2 + 2D cos 9 



df 



(15) 



In the above equation, the magnetization components are given by (m x ,m y ,m z ) = 
J (sin & cos 0', sin & sin 0', cos 9')f{9', 0', t) sin 0Wd0'. 

Now, consider an initial state prepared by sampling independently for each of the N 
spins the angle uniformly over [0, 27r] and the angle 9 uniformly over an arbitrary interval 
symmetric about tt/2. Such a state will have the distribution 

m<j>,Q) = ±p(0), (16) 
with p(9), the distribution for 9, given by 

if 9 e [f - a, f + a] , 



2 sin a 



P (0) 



(17) 







otherwise. 



We call such an initial state a waterbag state, in analogy with a similar form of an initial 
state studied in the context of the HMF model. It is easily verified that this non-magnetic 
state has the energy 
D 



sin 2 a, 



and that the state is stationary under the Vlasov dynamics (1151) . Let us proceed to analyze 
the dynamical stability of such a state. Such stability analysis in the context of the HMF 
model was pursued in [TB]. Here, we closely follow the treatment adopted in [23]. We 
consider finite but large system size N, and linearize the Vlasov equation (ITS"]) with respect 
to finite-size fluctuations 5f(9,<f),t) by expanding f(9,<j),t) as 



f(9, ( f>,t)^^p(9) + X5f(9,cf>,t). 



(19) 



Since the initial angles of the N spins are sampled independently, the small parameter A is 
of order 1/y/N. After linearization, Eq. ( ITS]) yields the following equation for 5f(0,(J),t): 



051 
dt 



m y cos i 



m r sin i 



1 dp{9) 



2D cos 9 



35 f 



(20) 



2tt d9 

where now m x and m y are linear in 5f. 

To study the linear dynamics, we note that at long times, it is dominated by the mode 
corresponding to the largest eigenfrequency uj of the linearized equation (1201 . Since the 
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perturbation Sf(9,4>,t) is 27r-periodic in </>, one has the following expansion in terms of the 
Fourier modes gk(9,oj), which is valid after a short initial transient (See [7J, Chapter 6): 

6f(0,<l>,t) = Y l 9k{0,u>)S kM >. (21) 



A' 

The magnetization components are given by 



771a; — 7T J sin 2 9' d9' (g^(9', to) + gi (9',uj))e luJt , (22) 

m y = -m J sm 2 9'd9'{g-i{9',cu)- gi {9',u))e iu3t . (23) 

It thus follows that the relevant eigenmodes of Eq. (1201) have k = ±1. Using Eq. (I2ip . 
and the above expressions for m x and m y in Eq. (TSUI) , we find that the coefficients #±1(0,0;) 
satisfy 

n/2+a 

9±i = dP }P ^r, 1 -irr—^ I #±i(0>)sin 2 0W. (24) 



d9 2(2Dcos9±u) 

ir/2-a 

On multiplying both sides of the above equation by sin 2 9 and then integrating over 9, we 
get 



/±(l-#±) = 0, (25) 



where 



and 



7r/2+a 

#±i(0,o;) sin 2 0dfl, (26) 

7r/2-a 



M6) S[n2() d9 (97) 
K± ~ U- a -dF2(2Dcos9±u) de - (27) 

Since I± 7^ 0, it then follows from Eq. (|25|) that the frequency u is given by the condition 

K± = l. (28) 
From Eq. ( 1T71) . we get 



f29l 

d# 2sinaL"V 2 ' J ~ \2 ' VJ' v ; 

which, together with Eqs. ( 1271) and ( 1281) . give the following expression for the largest 
eigenfrequency u: 

uj 2 = AD 2 sin 2 a — D cos 2 a. (30) 

Expressing the parameter a in the above equation in terms of the energy e in Eq. (fTBl) finally 
yields 

o; 2 = e(3 + 12D) — D. (31) 
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Figure 1. Dynamical phase diagram of our model in the (e,D) plane for the choice of the 
waterbag initial state, Eq. (fT6|) . The line to the right is the thermodynamic phase boundary 
D(e c ) between the magnetic and the non-magnetic phase, and is given by Eq. (|13j) with 
P c obtained by solving Eq. (fT2|) . The line to the left gives the dynamical phase boundary 
D(e*), and is given by Eq. ([32]). 
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We thus see that the frequency u is real for e > e*, given by 
D 



- 3 + 12/r (32) 

see Fig. (P). Therefore, unstable modes do not exist in this energy range so that the waterbag 
state ffl6|) is linearly stable. Hence, QSS is observed. In this case, in a finite system, such a 
state eventually relaxes to Boltzmann-Gibbs equilibrium on a timescale over which non-linear 
correction terms should be added to the Vlasov equation [3] . 

On the other hand, for e < e*, the waterbag state is unstable. Consequently, the 
perturbation Sf(9,4>,t) grows exponentially fast towards Boltzmann-Gibbs equilibrium. 
Below e*, on setting cu 2 = —O, 2 with real f2, we get 5f(9,(p,t) = Ae ±l( ^ +nt , where A is a 
constant. Consequently, the average magnetization behaves as 

m(t) ^e m ; e < e* , (33) 



before it relaxes to the equilibrium value. From the above equation, it follows that for e < e*, 
the relaxation timescale t(N) over which the magnetization acquires the equilibrium value 
of Ofl) scales as IniV. 



S 
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Figure 2. (a) Time evolution of average magnetization m(t) in the unstable phase 
with energy density e = 0.0603, the parameter D = 15, and for systems of size N = 
10000, 100000, 500000 (left to right), (b) Data collapse for the scaled magnetization VNm(t) 
as a function of £, in accordance with Eq. (|33[) . The black line represents the function 
s(t) ~ e nt , with f2 obtained from Eq. ([3T]) by substituting uj 2 = — £1 2 . Data averaging varies 
between 2 x 10 4 histories for the smallest system and 500 histories for the largest one. 

5. Numerical simulations 

In order to verify these features, we performed numerical simulations of the dynamics by 
integrating the equations of motion (@|, ((5]), fl6]) by using a fourth-order Runge Kutta method 
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tN- 1 - 7 

Figure 3. Average magnetization m{t) as a function of tN~ 1J in the stable phase 
with energy density e = 0.24, the parameter D — 15, and for systems of size N = 
300, 1000,3000,5000 (top to bottom). Data averaging varies between 300 histories for the 
smallest system and 25 histories for the largest one. The figure suggests a QSS life-time 
t(N) r^N 1 - 7 . 



with time step equal to 0.01. For e < e*, the results presented in Fig. [2](a) show that the 
magnetization grows fast towards equilibrium. On scaling the magnetization by vTV, Fig. 
EJ^b) shows a very good scaling collapse in accordance with the exponential growth predicted 
by Eq. ( 1331) . The growth rate Q is in agreement with that obtained from Eq. ( 13TT) by 
substituting u 2 = —Q 2 . For energies e* < e < e c , when the waterbag state (fl~6l) is linearly 
stable, Fig. [3] suggests a much longer relaxation time t(N) ~ N 1 - 7 . A similar scaling of the 
QSS relaxation time was also observed in the HMF model |18j . 

6. Conclusions 

In conclusion, we addressed the ubiquity of non-Boltzmann quasistationary states (QSS) 
during relaxation of long-range systems. This was done by studying an anisotropic 
Heisenberg model of globally coupled classical spins evolving under classical spin dynamics. 
Quasistationary states have earlier been shown to occur in long-range interacting systems 
composed of particles (inertial rotors) which are evolving under particle dynamics dictated 
by the underlying Hamiltonian. Thus, our model provides a different possible setting for 
the occurrence of QSS under spin dynamics. By analyzing the Vlasov equation for the time 
evolution of the phase space distribution, we demonstrated that in this model, relaxation of 



10 



a class of initial states in certain energy interval proceeds through intermediate QSS. These 
states have a lifetime that grows algebraically with the system size. This further establishes 
the possibility of long-range systems to exhibit quasistationarity under a broader class of 
dynamical processes. 

7. Acknowledgments 

We thank A. Bar, O. Cohen, T. Dauxois, O. Hirschberg, S. Levit and S. Ruffo for helpful 
discussions and comments on the manuscript. The support of the Israel Science Foundation 
(ISF) and the Minerva Foundation with funding from the Federal German Ministry for 
Education and Research is gratefully acknowledged. 

References 

[1] Dauxois T, Ruffo S, Arimondo E and Wilkens M (ed) 2002 Dynamics and Thermodynamics of Systems 
with Long-Range Interactions (Lecture Notes in Physics vol. 602) (Springer- Vcrlag, Berlin). 

[2] Campa A, Giansanti A, Morigi G and Sylos Labini F (ed) 2008 Dynamics and Thermodynamics of 
Systems with Long-Range Interactions: Theory and Experiments (AIP Conference Proceedings vol. 
970) (Melville, New York). 

[3] Campa A, Dauxois T and Ruffo S 2009 Phys. Rep. 480 57. 

[4] Dauxois T, Ruffo S and Cugliandolo L F (ed) 2010 Long-Range Interacting Systems (Oxford University 

Press, New York). 
[5] Bouchet F, Gupta S and Mukamcl D 2010 Physica A 389 4389. 

[6] Mukamcl D 2010 Long-Range Interacting Systems ed T Dauxois, S Ruffo and L F Cugliandolo (Oxford 

University Press, New York) 33. 
[7] Nicholson D R 1992 Introduction to Plasma Physics (Krieger Publishing Company, Florida). 
[8] Landau L D and Lifshitz E M 1960 Electrodynamics of Continuous Media (Pergamon, London). 
[9] Padmanabhan T 1990 Phys. Rep. 188 285. 
[10] Chavanis P H 2002 Dynamics and Thermodynamics of Systems with Long-range Interactions (Lecture 

Notes in Physics vol. 602) ed T Dauxois, S Ruffo, E Arimondo and M Wilkens (Springer- Vcrlag, 

Berlin) . 

[11] Lynden-Bcll D and Wood R 1968 Mon. Not. R. Astron. Soc. 138 495. 
[12] Thirring W 1970 Z. Phys. 235 339. 

[13] Barre J, Mukamel D and Ruffo S 2001 Phys. Rev. Lett. 87 030601. 
[14] Mukamel D, Ruffo S and Schreiber N 2005 Phys. Rev. Lett. 95 240604. 
[15] Bouchet F and Barre J 2005 J. Stat. Phys. 118 1073. 

[16] Bouchet F, Dauxois T, Mukamel D and Ruffo S 2008 Phys. Rev. E 77 011125. 
[17] Antoni M and Ruffo S 1995 Phys. Rev. E 52 2361. 

[18] Yamaguchi Y Y, Barre J, Bouchet F, Dauxois T and Ruffo S 2004 Physica A 337 36. 

[19] Campa A, Giansanti A and Morelli G 2007 Phys. Rev. E 76 041117. 

[20] Joyce M and Worrakitpoonpon T 2010 J. Stat. Mech.: Theory Exp. P10012. 

[21] Teles T N, Levin Y, Pakter R and Rizzato F B 2010 J. Stat. Mech.: Theory Exp. P05007. 

[22] Gabrielli A, Joyce M and Marcos B 2010 Phys. Rev. Lett. 105 210602. ' 

[23] Jain K, Bouchet F and Mukamel D 2007 J. Stat. Mech.: Theory Exp. P11008. 

[24] Latora V, Rapisarda A and Ruffo S 1999 Phys. Rev. Lett. 83 2104. 



11 



[25] Latora V, Rapisarda A and Tsallis C 2001 Phys. Rev. E 64 056134. 

[26] Bouchet F and Dauxois T 2005 Phys. Rev. E 72 045103(R). 

[27] Nobre F D and Tsallis C 2003 Phys. Rev. E 68 036115. 

[28] Baldovin F and Orlandini E 2006 Phys. Rev. Lett. 96 240602. 

[29] Baldovin F and Orlandini E 2006 Phys. Rev. Lett. 97 100601. 

[30] Baldovin F, Chavanis P H and Orlandini E 2009 Phys. Rev. E 79 011102. 

[31] Gupta S and Mukamel D 2010 Phys. Rev. Lett. 105 040602. 

[32] Gupta S and Mukamel D 2010 J. Stat. Mech.: Theory Exp. P08026. 

[33] Chavanis P H, Baldovin F and Orlandini E 2010 \arXiv:1009.5603\ 

[34] Kac M, Uhlenbeck G E and Hemmer P C 1963 J. Math. Phys. 4 216. 

[35] Mermin N D 1967 J. Math. Phys. 8 1061. 

[36] Barre J, Mukamel D and Ruffo S 2002 Dynamics and Thermodynamics of Systems with Long-range 
Interactions (Lecture Notes in Physics vol. 602) ed T Dauxois, S Ruffo, E Arimondo and M Wilkens 
(Springer- Verlag, Berlin) . 

[37] Braun W and Hepp K 1977 Comm. Math. Phys. 56 101. 



12 



